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1.  Introduction 

A  dominant  process  for  radiation  belt  electrons  is  pitch 
angle  diffusion,  caused  by  cyclotron  and  Landau  resonance 
with  a  variety  of  plasma  waves.  The  resulting  scattering 
into  the  loss  cone  and  precipitation  into  the  atmosphere  is 
a  major  loss  mechanism,  especially  inside  of  geosynchro¬ 
nous  orbit  Neglecting  all  other  processes,  the  relaxation 
of  a  population  injected  with  an  arbitrary  pitch  angle 
distribution  can  be  considered  to  occur  on  a  sequence 
of  time  scales  given  by  the  ordered  eigenvalues  of  the 
bounce-averaged  pitch  angle  diffusion  operator  (eg., 
Schulz  and  Lanzerotti,  1974;  see  also,  O’Brien  et  al., 
2008).  Although  the  long  time  behavior  is  described  by 
the  single,  longest  timescale,  the  time  required  for  the 
higher  eigenmodes  to  decay  can  be  significant,  and  can 
lead  to  complications  in  interpreting  short  sequences 
of  observations  (Kress  et  aL,  2007,  Baker  et  al.,  2007).  In 
this  paper,  we  generally  refer  to  the  slowest  diffusion  time 
scale  as  "the  lifetime.” 

Strictly,  this  eigenvalue  problem  depends  on  the  values 
of  the  diffusion  coefficients  at  all  values  of  equatorial  pitch 
angle  a0.  These  coefficients  can  vary  widely  with  a0.  as 
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different  resonances  come  into  play  depending  on  the 
local  values  of  particle  pitch  angle,  wave  refractive  index, 
and  geomagnetic  field.  Typically,  for  a  given  wave  mode, 
diffusion  by  cyclotron  resonances  is  the  largest  for  small 
oco  while  diffusion  by  the  Landau  resonance  is  the  largest 
for  large  a0,  leading  to  a  pronounced  minimum  at  some 
intermediate  <x0  (although  this  may  be  masked  by 
combining  different  wave  modes  with  diffusion  minima 
located  at  different  a0  values).  Lyons  and  Thorne  (1973) 
gave  an  integral  equation  for  simultaneously  determining 
the  lowest  eigenvalue  and  eigenfunction,  or  the  lifetime 
and  relaxed  pitch  angle  distribution,  for  an  arbitrary 
profile  of  D(ao).  Albert  (1994)  reformulated  this  as  a 
boundary  value  problem  for  a  set  of  ordinary  differential 
equations,  readily  solved  with  "shooting"  method  (Press 
et  al.,  1992).  Both  approaches  will  be  used  below. 

Shprits  et  al.  (2006)  proposed  a  simple,  convenient 
method  of  estimating  the  global  diffusion  lifetime,  namely 
as  the  inverse  of  the  diffusion  coefficient  evaluated  at  the 
edge  of  the  loss  cone  (denoted  by  Dt).  This  follows  in 
the  footsteps  of  several  earlier  papers  which  related  the 
lowest  eigenvalue  to  DL  for  simple  model  problems  (e.g., 
Kennel  and  Petschek,  1966;  Roberts,  1969;  Schulz 
and  Lanzerotti,  1974).  Shprits  et  al.  (2006)  argued  that 
away  from  the  loss  cone,  small  diffusion  coefficients 
are  compensated  by  the  development  of  large  gradients 
in  the  particle  distribution.  A  series  of  numerical  experi¬ 
ments  with  model  diffusion  coefficients  showed  reasonable 
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agreement  for  D(a0)  profiles  with  moderate  internal 
minima.  They  also  noted  that  large  minima  lead  to  slow 
relaxation  of  the  pitch  angle  profile,  perhaps  slower  than 
of  physical  interest,  while  the  behavior  on  shorter  time 
scales  did  seem  to  be  controlled  by  the  loss  cone  diffusion 
rate.  They  concluded  that  internal  minima  in  D  were 
subordinate  to  the  loss  cone  value  as  long  as  these 
diffusion  coefficients  were  within  a  factor  of  10. 

On  the  other  hand,  one  might  expect  the  decay  rate  to 
be  dictated  by  the  minimum  value  of  the  diffusion  rate, 
termed  a  "bottleneck”  by  Schulz  (1991)  who  tied  them  to 
"shoulders”  in  the  distribution  function  /(a 0).  Lichtenberg 
and  Lieberman  (1983),  discussing  transport  with  variable 
diffusion  coefficients  in  a  different  context,  also  noted  that 
the  gradients  becomes  steeper  where  the  diffusion  is 
slower  but  concluded  from  the  continuity  of  flow  that  the 
effective  diffusion  rate  is  given  by  the  harmonic  average, 
((1  /D))~} ,  and  that  this  is  dominated  by  regions  where  D  is 
the  smallest.  This  seems  consistent  with  the  numerical 
results  of  Abel  and  Thorne  (1998),  who  found  that 
diffusion  due  to  lightning  whistlers  and  VLF  transmitters, 
while  negligible  compared  to  Coulomb  collisions  and 
hiss  rates  at  the  loss  cone,  greatly  reduce  the  lifetime  of 
500  keV  electrons  at  L  —  2.4  by  smoothing  the  deep 
minimum  that  would  otherwise  occur  at  a0  %  65°. 
Similarly  a  numerical  study  of  diffusion  due  to  combined 
whistler  mode  hiss  and  EM1C  waves  (Albert,  2003)  found 
the  electron  lifetimes  to  be  strongly  affected  by  the  hiss 
amplitudes  which  determined  the  value  (and  location)  of 
the  minimum  in  D,  even  though  hiss  was  negligible  at 
theloss  cone  itself. 

To  explore  and  reconcile  these  two  viewpoints,  we 
start  with  exact  expressions  for  the  decay  rates  and  derive 
an  integral  approximation  which  reflects  aspects  of  both 
We  also  evaluate  a  model  that  admits  a  (nearly)  exact 
analytical  treatment  Section  3  tests  the  approximation 
and  several  simpler  expressions  against  an  analytically 
solvable  model,  and  Section  4  does  the  same  with  a  series 
of  numerical  experiments  paralleling  those  of  Shprits  et  al. 
(2006).  We  conclude  that  the  integral  approximation  is  an 
accurate  and  reliable  estimate  for  a  wide  variety  of 
diffusion  models. 


2.  General  expressions  for  lifetimes 


We  start  with  the  one-dimensional  bounce-averaged 
pitch  angle  diffusion  equation, 


0/  1  a 

6t  T  sin  2a0  03(o 


DT  sin  2a0 


+  S, 


(1) 


where  ao  is  the  equatorial  pitch  angle,  D  is  the  cor¬ 
responding  pitch  angle  diffusion  coefficient,  S  is  a 
stationary  source  term,  and  T  ^  1  30  —  0.56  sin  ao  is  the 
normalized  bounce  time.  The  boundary  conditions  are /  = 
0  at  the  loss  cone  angle  cll  and  df  /dcto  =  0  at  a0  =  n/2. 

Following  Lyons  and  Thorne  (1973),  we  write  /(a0,t)  = 
/o(2o)  +  g(o(o)exp(-t/T),  where /0  solves  Eq.  (1)  in  steady 
state,  to  get 


-gT  sin2a0 

T 


d 

da0 


^DT  sin2a0 


(2) 


Integrating  from  the  loss  cone  angle  aL  to  n/2  and  taking 
the  boundary  conditions  into  account  gives 


1  fn/2  dj? 

-/  gTsin2a'da'  =  D/.T/.  sin  2  aL:r^- 
T  J  Xl  doto 


(3) 


where  DL  =  D(otL )  and  TL  =  T(aL).  Eq.  (3)  is  the  same  as 
Eq.  (13)  of  Lyons  and  Thorne.  More  generally,  integrating 
from  an  arbitrary  value  of  a  to  7t/2  gives 


1  fn/2  dg 

-j  gT  sin  2a'  da'  =  DT  sin  2a  — (4) 

Dividing  Eq.  (3)  by  Eq.  (4),  solving  for  dg/da,  and 
integrating  from  ctL  to  a0  gives 


d 

DTsin2a 


where 


(5) 


A(a)  =  1 


/jigTsin2a'da' 
I*!2  &T  sin  2a' da' 


(6) 


Eq.  (5)  is  the  same  as  Eq.  (12)  of  Lyons  and  Thorne  (1973). 
Combining  this  with  Eq.  (3)  gives 


rn/2  ^ c 

=  /  daoTsin2o(o  / 

Jy-i  J*l 


da 


DT  sin  2a 


/1(a). 


(7) 


We  denote  by  t0  the  value  of  r  obtained  with  the 
approximation  T(a0)=l.  Making  this  approximation  in 
Eq.  (7)  and  reversing  the  order  of  integration  according  to 

51 2  doc0  C  da  =  51 2  da  /* '2  da0  yields 


:  t0 


-f 

JrxL 


71/2  ,  cos  a 

ds-— - — 
2D  sin  a 


/1(a). 


(8) 


3.  An  estimate 


We  now  consider  the  situation  where  D  has  a  deep 
minimum  at  some  value  a  =  a*.  As  a  decreases  from  n/ 2, 
the  left-hand  side  of  Eq.  (4)  increases.  If  D  becomes  small 
at  a  =  a*,  then  dg/da  must  become  large  there,  so  thatg  is 
much  larger  for  a>a,  than  for  a<a*.  If  g(a)sin2a  is 
approximated  as  0  for  a  <a,  and  constant  for  a  >  a„  Eq.  (6) 
gives 

{1,  a^a*, 

a>a.  (9) 

so  that 


T 


da 


cos  a 
2D  sin  a 


da 


cos  a  a  —  a* 

2D  sin  a  (7i/2)  —  a^  * 


(10) 


Because  the  second  integrand  is  zero  at  a*,  where 
l/(Dsina)  is  assumed  to  have  a  steep  maximum,  the 
second  integral  will  be  less  than  the  first  integral,  finally 
yielding  the  estimate 


da 


cos  a 
2D  sin  a ' 


(ID 


This  formula  clearly  suggests  that  t  will  be  determined 
mostly  by  the  minimum  value  of  2Dtana0. 
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This  partially  vindicates  both  simple  estimates  dis¬ 
cussed  above.  If  the  minimum  of  D  is  deep  enough,  it  will 
largely  determine  the  minimum  of  Dtana0.  If  D(xo)  is 
relatively  flat,  the  minimum  of  Dtana0  will  occur  at  the 
loss  cone.  If  the  minimum  of  D  occurs  near  the  loss  cone, 
the  two  concepts  coincide. 

4.  An  analytically  solvable  model 


For  some  simple  analytical  forms  of  D(xo),  and 
neglecting  the  variation  of  T(x0),  it  is  possible  to  write 
the  complete  set  of  eigenvalues  and  eigenfunctions  of  the 
diffusion  operator  in  closed  form.  As  given  by  Schulz 
and  Lanzerotti  (1974),  with  slight  change  of  notation, 
the  model  D  =  D0  cos2a  x0/ sin2  x0  with  a<  1  yields  (as  we 
have  verified)  the  eigenvalues  and  (unnormalized)  eigen¬ 
functions 


(i  -  arn-in 

An  —  -» 


g„(a o>  = 


COS  z~zoXL 
COSa  Oti 

cos'7  x 


Jv 


n  —  0, 1, 2, . . . 

cos1_CT  x 
cos1-*7  aj* 


(12) 


where  v  =  <x/(  1  —  rr)  and  Kvn  is  the  nth  zero  of  the  Bessel 
function  Jv.  Remarkably,  the  expressions  remain  valid  for 
(7<0,  in  which  case  D  has  an  internal  minimum  at 
sin2  a0  =  1/(1  -  a);  the  minimum  value  is  D0(  1  —  cr) 
(-rr/(l  -  (T)f  Also,  rr<0  guarantees  that  the  eigenfunc¬ 
tions  obey  the  boundary  conditions  associated  with 
Eq  (1). 

With  o  =  -1  Jv(z)  =  ( 2/7tz )1/2  cosz  (Gradshteyn  and 
Ryzhik,  1980),  and  the  eigenvalues  and  eigenfunctions 
become 


7t2(  1  +  2  n)2 
cos4  (XL 


Do, 


n  =  0,1,2,. 


gn(a0)  -  cos 


7t(  1  +  2n)cos2  x0 
2  cos2  xL 


(13) 


These  results  agree  closely  with  full  numerical  solutions 
based  on  Albert  (1994),  and  essentially  perfectly  with 
numerical  solutions  setting  T(a0)  =  1 
The  inverse  of  the  lowest  eigenvalue  is 


cos4  xL 

T,)=^w 

From  Eq.  (11),  the  estimate  x.  is 


(14) 


D  O  1/t0 


- g  a  1  /t. 

9  (T=1)  +  Di 

•  1 A  X 


Fig.  1.  Diffusion  coefficient  D  and  pitch  angle  profile  g  vs.  equatorial 
pitch  angle  oc0  for  the  analytical  model  discussed  in  the  text  Also  shown 
are  the  decay  lifetime  r  and  several  estimates. 


loss  cone, 

1  cos2  aL  sin11  at 

w=  w  •  (18) 

differs  by  the  factor  ~10tan2  a l,  which  can  be  significantly 
less  than  1. 

Fig.  1  shows  D(ao)  and  the  shape  of  g(a0),  calculated 
with  Xl  corresponding  to  L  =  4  in  a  dipole  magnetic  field. 
The  approximation  T(a0)  =  1  has  virtually  no  effect  on  the 
shape  of  g.  Also  shown  are  the  value  of  x  from  a  full 
computation,  the  analytic  value  x0,  and  the  estimate  x,,  all 
of  which  agree  quite  closely,  as  well  as  the  simpler 
estimates  Di  and  Dmin.  Both  simpler  estimates  are  seen  to 
give  significantly  poorer  agreement  than  x.. 


5.  Numerical  experiments 


cos4  XL 
~SDq~ 


(15) 


which  agrees  extremely  well  with  x0,  even  though  the 
minimum  of  this  model  D{x0)  is  not  sharp  enough  to  cause 
a  pronounced  "shoulder”  in  g(a0)  as  assumed.  The  inverse 
minimum  values 


1  1 

^min  4Do 
and 


(16) 


1  3/3  ^  1 

(2Dtana0)min  ”  32 D0  6D0 


(17) 


also  usually  agree  with  x0  to  within  about  a  factor  of  2  for 
plausible  values  of  xL.  The  inverse  value  at  the  edge  of  the 


To  further  test  the  estimate  given  by  Eq.  (11 ),  we  repeat 
the  numerical  experiments  of  Shprits  et  al  (2006).  These 
were  based  on  pitch  angle  diffusion  coefficients  for  1  MeV 
electrons  at  L  =  4  due  to  lower-band  chorus  waves, 
calculated  using  the  parallel-propagation  formulation 
of  Summers  (2005).  (The  parameters  used  were  wM  = 
0  35£>e,  Sto  =  0  15f2e,  c/;Lc  =  0  05f2e,  ojuc  =  0  65f2ff,  and 
wpe/Qe  =  4A5.)  The  D  values  were  combined  with 
artificial  step-function  decreases  at  various  a0  ranges. 
For  each  modified  version  of  D,  we  computed  the  lowest 
eigenmode  g(a0)  and  lifetime  x  according  to  the  shooting 
code  algorithm  of  Albert  (1994),  and  repeated  the 
calculation  with  the  approximation  T(a0)  =  1  to  obtain 
x0.  The  lifetimes  were  compared  to  the  approximation  x* 
and  the  simple  estimates  D^1  and  D~]n. 
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Fig.  2.  Diffusion  coefficients  based  on  parallel-propagating  chorus  waves  combined  with  step  decreases  at  various  locations  Also  shown  are  the  lifetime 
and  various  estimates;  the  symbols  have  the  same  meaning  as  in  Fig,  1. 
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Fig.  3.  Same  as  Fig.  2  but  with  step  decreases  of  various  depths. 


Fig.  2  corresponds  to  Fig.  3  of  Shprits  et  al.  (2006), 
where  a  step  of  depth  DL/ 10  and  width  20°  was  imposed 
at  large,  intermediate,  and  small  values  of  20-  In  all  cases, 
both  r  and  g(a0)  agree  with  the  results  found  by  Shprits 
et  al.  by  directly  evolving  the  diffusion  equation  in 
time  numerically.  The  estimate  Dj~J  works  well  if  the 
step  is  at  large  or  small  a0.  or  absent,  but  starts  to 
break  down  when  the  step  is  at  intermediate  oc0. 


The  estimate  D“]n  is  generally  poor,  but  t*  is  consistently 
good.  Actually,  t.  is  a  better  estimate  of  To  than 
of  t,  since  both  assume  T(a0)  =  1,  but  the  difference  is 
small. 

Fig.  3,  corresponding  to  Fig.  4  of  Shprits  et  al.  (2006), 
shows  a  similar  comparison  with  the  step  decrease  held  at 
intermediate  a0  but  with  varying  depth  The  value  1  is  a 
good  approximation  of  r  for  small  and  moderate  step 
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Fig.  4.  Same  as  Fig.  2  but  with  step  decreases  of  various  widths. 
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Fig-  5.  Same  as  Fig.  2  but  with  large  step  decreases  at  large  and  small  a0. 


decreases,  but  breaks  down  for  large  steps.  D~|n  is  a  poor 
estimate  of  x,  but  r„  is  reliably  good. 

In  Fig.  4,  the  width  of  the  step  is  varied.  Once  again,  t 
and  T0  are  approximated  poorly  by  Dm jn  and  well  by  DL 1 
for  narrow  and  moderate  steps,  but  the  agreement  breaks 
down  for  wide  steps,  while  x,  provides  a  uniformly  good 
estimate. 

Finally,  Fig.  5  shows  much  deeper  steps  (DL/100  and 
Dt/1000).  When  the  step  is  at  small  a0,  all  the  values  are 


close  together.  However,  when  the  step  is  at  large  a0,  only  r. 
is  a  good  approximation  to  t  and  t0. 

6.  Conclusions 

We  have  derived  a  new  estimate  of  the  lowest  pitch 
angle  diffusion  timescale  x,  namely  x,  of  Eq.  (11),  based 
on  the  integral  of  1/Dtana,  and  verified  its  accuracy 
with  an  analytically  solvable  model  and  with  a  variety  of 
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numerical  experiments.  The  approximation  is  simpler  and 
more  physically  transparent  than  the  full  calculation,  and 
allows  convenient  estimates  of  changing  various  wave 
parameters  (such  as  amplitude  or  latitudinal  extent)  by 
considering  D  at  only  a  certain  values  of  a0-  The  exact 
expression  z  given  by  Eq  (7).  however,  is  not  difficult  to 
compute  (as  in  Albert.  1994),  and  is  the  most  rigorous. 

Based  on  their  numerical  experiments,  Shprits  et  al. 
(2006)  concluded  that  the  estimate  of  z  by  D[]  was  often 
good  but  broke  down  for  step  function  minima  in  D  less 
than  D/./10.  This  makes  sense  in  light  of  Eq.  (11),  since 
the  integral  will  be  dominated  by  the  minimum  value 
of  D  tan  a,  and  that  minimum  will  only  be  driven  by  a 
minimum  in  D  itself  where  (for  moderate  a),  D<Dt  tan  aL. 
Since  typically  tan  (Xl  %  1/10,  a  minimum  in  D  would  have 
to  be  less  than  0.1  Di  to  be  dominant  A  similar  conclusion 
can  be  drawn  from  the  analytical  results  of  Section  4: 
there,  z  becomes  substantially  larger  than  D^1  when 
7t sin  ai<l.  which  implies  Dmin<(4/n2)DL^  QADL.  As 
noted  by  Shprits  et  al.  (2007),  actual  diffusion  coefficients 
can  exhibit  deep  minima,  such  as  when  EM1C  waves  are 
present  in  combination  with  whistler  mode  chorus  waves. 

It  is  also  true  that  the  resulting  lowest  order  lifetimes 
may  be  longer  than  the  duration  of  the  waves,  so  that  the 
particle  population  may  not  have  time  to  settle  into 
the  lowest  eigenmode.  In  these  situations,  and  especially 
considering  energy  diffusion  and  radial  transport,  no 
single  timescale  may  be  sufficient  to  describe  the 
dynamics  in  a  meaningful  way.  Fortunately,  the  develop¬ 
ment  of  2D  and  3D  codes  to  advance  the  time-dependent 
diffusion  equation  is  well  under  way. 
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